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In this paper, we present the consistent evolution of short-period exoplanets coupling the tidal and gravothermal evo- 
lution of the planet. Contrarily to previous similar studies, our calculations are based on the complete tidal evolution 
equations of the Hut (1981) model, valid at any order in eccentricity, obliquity and spin. We demonstrate, both analyt- 
ically and numerically, that, except if the system was formed with a nearly circular orbit (e;^0.2), solving consistently 
the complete tidal equations is mandatory to derive correct tidal evolution histories. We show that calculations based 
on tidal models truncated at second order in eccentricity, as done in all previous studies, lead to quantitatively, and 
sometimes even qualitatively, erroneous tidal evolutions. As a consequence, tidal energy dissipation rates are severely 
underestimated in all these calculations and the characteristic timescales for the various orbital parameters evolutions 
can be wrong by up to three orders in magnitude. Such discrepancies can by no means be justified by invoking the 
uncertainty in the tidal quality factors. 

Based on these complete, consistent calculations, we revisit the viability of the tidal heating hypothesis to explain 
the anomalously large radius of transiting giant planets. We show that, even though tidal dissipation does provide a 
substantial contribution to the planet's heat budget and can explain some of the moderately bloated hot-Jupiters, this 
mechanism can not explain alone the properties of the most inflated objects, including HD 209 458 b. Indeed, solving the 
complete tidal equations shows that enhanced tidal dissipation and thus orbit circularization occur too early during the 
planet's evolution to provide enough extra energy at the present epoch. In that case, either a third, so far undetected, 
low-mass companion must be present to keep exciting the eccentricity of the giant planet, or other mechanisms, such as 
stellar irradiation induced surface winds dissipating in the planet's tidal bulges and thus reaching the convective layers, 
or inefficient flux transport by convection in the planet's interior must be invoked, together with tidal dissipation, to 
provide all the pieces of the abnormally large exoplanet puzzle. 
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1. Introduction 

Gravitational tides have marked out the history of science 
and astrophysics since the first assessment by Seleucus of 
Seleucia of the relation between the height of the tides and 
the position of the moon and the Sun in the second cen- 
• • tury BC. Modern astrophysics extended the study of grav- 
. !^ itational tides in an impressive variety of contexts from the 
k> synchronization of the moon and other satellites to the evo- 
lution of close binary stars and even the disruption of galax- 
d ies. 

The recent discoveries of short period extrasolar plan- 
etary systems and the determination of the anomalously 
large radius of some giant close-in exoplanets revived the 
need for a theory of planetary tides covering a wider va- 
riety of orbital configurations than previously encountered 
for the case of our own solar system planets. In particu- 
lar, the orbital evolution of planetary systems such as HD 
80 606, with an orbital eccentricity of 0.9337 (Naef et al. 



tidal models limited to the case of zero or vanishing eccen- 
tricity and obl i quity such as in the m ode ls of e.g.|Goldreich 
'fc Soter| ( |196 6[), [Jackson et al.|(|2008[) and|Ferraz-Mello et al. 
(2008r 



Following Bodenheimer et al. (2001) and Gu et al. 



(2003), attempts have been made to explain the observed 
large radius of some transiting close-in gas giant exoplanets 
the so-called "Hot Jupiters" - by means of tidal hea ting 



20011, and XO-3 with a stellar obhquity > 37.3 ± 3.7 
deg (Winn et al. 2009b I, cannot be properly treated with 
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( jJackson et al.|2008[plller et al.|2069{[rbgui et al.|2009[ ). All 

these models, however, use tidal models truncated to low 
(second) order in eccentricity, in spite of initial eccentric- 
ities, as determined from the tidal evolution calculations, 
which can be as large as e = 0.8! According to these calcu- 
lations, a large eccentricity can remain long enough to lead 
to tidal energy dissipation in the planet's gaseous envelop 
(assuming a proper dissipation mechanism is at play in the 
deep convective layers) at a late epoch and then can explain 
the actual bloated radius of some observed planets. 

In the present paper, we revisit the viability of this tidal 
heatin g hypothesis, using an extended version of the [mit] 
( 1981 ) tidal evolution model, solving consistently the com- 
plete tidal equations, to any order in eccentricity and obliq- 
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uity, and coupling these latter with the gravothermal evo- 
lution of the irradiated planet. As will be shown in the 
paper, properly taking into account the full nature of the 
tidal equations severely modifies the planet's tidal and ther- 
mal evolution, compared with the aforementioned trun- 
cated calculations, leading to significantly different tidal 
heat rates and thus planet contraction rates. 

After introducing our model in f|2| we examine in de- 
tail in f|3] the relation between the constant time lag (At) 
in Hut's (and thus our) model and the usual tidal quality 
factor (Q) widely used in the literature. Constraints on A< 
from the study of the galilean satellites are also derived. In 
^ we demonstrate, with analytical arguments, that trun- 
cating the tidal equations at 2"^^ order in eccentricity leads 
to wrong tidal evolution histories, with sequences drasti- 
cally differing from the ones obtained when solving the com- 
plete equations. In Sjs] we compare our full thermal/orbital 
evolution calculations with similar studies based on a trun- 
cated and constant Q tidal model. These numerical com- 
parisons confirm and quantify the conclusions reached in 
^ namely that low order eccentricity models substantially 
underestimate the tidal evolution timescales for initially 
eccentric systems and thus lead to incorrect tidal energy 
contributions to the planet's energy balance. For instance, 
we show that tidal heating can not explain the radius of 
HD 209 458 b, for the present values of their orbital param- 
eters, contrarily to what has been claimed in pre vious calcu- 
lation s based on truncated eccentricity models ( Ibgui et al. 
2009 1. Finally, in [6j we apply our model to the case of some 
of the discovered bloated planets. We show that, although 
tidal heating can explain the presently observed radius of 
some moderately bloated hot Jupiters, as indeed suggested 
in some previous studies, tidal heating alone cannot explain 
all the anomalously large radii. Indeed, in these cases, ec- 
centricity damping occurs too early in the system's tidal 
evolution (assuming a genuine two-body planetary system) 
to lead to the present state of the planet's contraction. 



2. Model Description 

2.1. Internal evolution 

The main physics inputs (equations of state, internal com- 
position, irradiated atmosphere models, boundary condi- 
tions) used in the present calculations have been described 
in details in previous pa pers devoted to the evo l ution of 
extrasolar giant planets (liaraffe et al. 12003, 2008 



Leconte 

et al. 20091 and are only briefly outlined below. The evo- 
lution of the planet is based on a consistent treatment be- 
tween the outer non-grey irradiated atmospheric structure 
and the inner structure. The interior is composed primarily 
of a gaseous H/He envelope whose thermodynamic proper- 
ties are described by the Saumon-Chab rier-VanHorn equa- 
tion of state (EOS, Saumon et al.|1995 ) with a solar or non- 
solar enrichme nt in heavy elements described by the appro- 
priate EOS's (Baraffe et al. 2008). In the present calcula- 
tions, our fiducial model consists of a planet wi th a central 
core made up o f water, with the ANEOS EOS (Thompson 
h Lauson|1972 1 . A detailed analysis of the effects of differ- 
ent EOS's, core compositions and heavy material reparti- 
tions within the planet can be found in Baraffe et al.^(2008| , 
as well as a compariso n with models from o ther groups, in 
particular the ones by Fortney et al. (2007). 



Transiting planets are by definition very close to their 
host star (a < O.IAU). In that case, the stellar irradiation 



levels (Barman et al. |2001|) and thus the 


planet's evolu- 


tion 


Ouiliot et ai.||1996l Baraffe et al.||2003 


Burrows et al. 


2003 


Chabrier et al. 2004). We use a gri( 


of irradiated 



et al. (2001), computed for different levels of stellar irradi- 
ation relevant to the present study. For planets with a finite 
orbital eccentricity, the mean stellar flux received is given 
by: 



where Ri, and Tcff .* are the stellar radius and effective tem- 
perature, a the Stefan-Boltzmann constant, a and e the 
orbital semi-major axis and eccentricity and / is a geomet- 
rical factor (substituting / by 1 in Eq.[l]gives the mean flux 
received at the substellar point). The atmos pheric models 
were c omputed with / = 1/2 as described by [Barman et al. 

pooTl). 



2.2. Tidal Model 

We consider the gravitational tides raised by both the host 
star and the planet on each other and follow the tr aditional 
"viscous" approach of the equilibrium tide theory (Darwin 



be calculated exactly (e.g. 


Hut 


l!)M 


Neron de Surgy & 


Laskar 1997( |Levrard et al. 


|200 


I\ ICorreia & Laskar 2010 



value of the eccentricity and obliquity) 
1 da 4 a 
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where G is the gravitational constant, Wp is the planet's ro- 
tation rate, Mp and Rp its mass and radius, Sp its obliquity 
(the angle between the equatorial and orbital planes), with 



cosffp, and n the orbital mean motion. 



Kp = -fc2,pAtp 



GAfp2 

Rn 



Rr. 



(5) 



where fc2,p and Aip are the potential Love number of de- 
gree 2 and the constant time lag for the planet. The stellar 
parameters correspond to the same definitions, by simply 
switching the p and * indices. Similarly, the secular evolu- 
tion of the eccentricity is given by 

Ide 

e di 



= 11- 



" GA'UMp 
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with 



and 



2^ 8^ 
2^5 



(l-e2) 



15p2 , 15 4 I 5_ 6 
(1 - e2)13/2 



(7) 



(8) 



The terms proportional to Kp (resp. if^,) are due to the 
tides raised by the star (planet) on the planet (star). 
Finally, the evolution of the rotational state of each object 
i (i being p or ★) is given by: 



dCiWi 
At 



n 



1 



2x;N{e) 



(9) 



while the evolution of the obliquity obeys the equation 



— = smei- 

at CiCjj n 



{xi-r]{)n{e) — -2N{e) 



(10) 



where Ci is the principal moment of inertia of the de- 
formable body under consideration, 77; is the ratio of ro- 
tational over orbital angular momentum 



Mp + CiWi 



a'^ny/Y 



and 



l + 3e2 



(1 -e2)9/2 



(11) 



Up to this point, no assumption has been made on the 
objects themselves. As a result, Eqs. (pl)-([lO| are fully sym- 
metric in p and * indices and can be used directly to model 
binary stars or a planet-satellite system. For typical HD 
209 458 b-like parameters and Ai- values comparable to that 
inferred for Jupiter (see ^ , the planetary spin evolves to a 
coplanar (cp — 0) state with the equilibrium rotation rate 
value (setting dujp/dt = 0) 



N{e) 2x 
il{e) 1 + x 



Nje) 



n, 



(12) 



within a time scale r = Cp n^/Kp ^ 10^ yr. Since the age 
of the known transiting systems ranges from a few Myr 
to Gyr, we can safely assume that the planet is now in 
this state of pseudo synchronization. The evolution of the 
orbital parameters of the planet-star systems is thus fully 
defined by the system of differential equations, to be solved 
consistently, given by Eqs. ^ and ^ for a and e, Eqs. ^ 
and ((Tol) for the stellar rotation and Sp — 0, ujp — uj^q- 



The rate of tidal dissipation with in the planet in this state 
of ps eudo synchronization reads ( |Hut||198l Levrard et al. 
^OOTl 



-Etidcs — 2i^p 



Na{e) 



{ujp = Wcq) 



(13) 



The dissipated heat is deposited over the whole planet's 
interior. 



One can see from Eq. (A.26) in Appendix A that 



Eq. ( 13 1 is a special case of energy dissipation for a body in 
pseudo synchronous rotation as expected for fiuid objects 



0, Wn 



:j). For a rocky planet, the external grav- 



itational potential created by its permanent quadrupole 



moment can cause its locking into synchronous rotation 
{ujp = n) and the dissipation rate reads in that case 



E, 



tides 



2Kp X 



/ 1 + 

iV,(e) - 2A^(e) Xp + fi(e) 



(wp = n) 
(14) 



This equation fully agrees with Eq. (30) of Wisdom] ( 2008 ) 
who calculated it for a homogeneous, incompressible with a 
radial displacement Love number /12 = 5A;2/3. Note that our 
derivation does not require such an hypothesis and all the 
uncertainties in the radial distribution of material and its 
physical properties (e.g., density, co mpr essibility, elasticity) 
are lumped into the k2 parameter ( jLevrard„2008) . 

3. Relationship between the time lag At and the 
quality factor (Q) 

The aforcdcscribcd tidal model, which leads to exact tidal 
evolution equations in the viscous approximation, implies 
a constant time lag At. Neither the tid al quality factor 
(Q), or its counterpart, th e phase lag (e) ([Goldreich 



1963 



Goldreich fc Soter 19661 enter the dynamical evolution 
equations. Instead, the model is characterized by the time 
lag between the maximum of the tidal potential and the 
tidal bulge in each body, Atj, consid ered to be con stant 
duri ng the ev olution. As shown e.g. by Darwin ( 1908 ) (see 
also Greenberg 2009), this model is equivalent to consid- 
ering a body whose rheology entails Q~^{a) « e{a) (x a, 
where a/27T is the frequency of the tidal forcing. The actual 
rheology of giant gaseous planets being poorly constrained, 
this arbitrary choice based on the visco-elastic model has 
the advantages of (i) not introducing any discontinuity for 
vanishing tidal frequencies, as is the case for synchronous 
rotation, and (ii) to allow for a complete calculation of the 
tidal effect without any assumption on the eccentricity for 
an ideal viscoelastic bod)^; 



Indeed, as shown by Greenberg (2009), the frequency 
dependence of the phase lag of a perfect viscoelastic oscil- 
lator is given by 



tan(e) 



(15) 



where t is a viscous damping timescale and ujq the natu- 
ral frequency of the oscillator. In an incompressible gaseous 
body, the restoring force acting against the tidal deforma- 
tion is the self gravity of the body. Thus luq can be estimated 

through the free-fall time as 2tt/u>o w Y^fi ~ ^0 min- 
utes for Jupiter mean density. For tidal periods of several 
days luq 3> cr, and for weakly viscous fluid, the phase lag 
reads 

e((T) w = aAt (16) 

T 

which is the frequency dependence corresponding to the 
constant time lag model. 

On the contrary, constant-Q models described by 



ierl 

McUo et al. (2008T were derived using perturbative devel 
opments of Kepler equations of motion both in eccentric- 
ity and inclination. Such Fourier decomposition is in fact 



Goldreich & Soter (1966), Jackson et al. 



( 


2008 


1, 


Ferraz- 
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necessary in an "lag and add" approach with a given fre- 
quency dependence of the phase lag (e((T)). Indeed, in this 
approach, one must first separate the forcing potential in 
terms with a defined frequency before lagging them with 



the c hosen e(a) (See Ferraz-Mello et al 



2008 



2009 1 . As a result they can only be used in the e ^ 1 and 



Ei <^ 1 limit. 

The time lag At can be linked to the reduced quality 
factor Q' = 3Q/2k2, chosen so that Q' — Q for an ho- 
mogeneous sphere (^2 = 3/2). Indeed, one must remember 
that the phase lag, e(cr), induced by the tidal dissipative ef- 
fects, is twice the geometrical lag angle, S{(t), between the 
maximum of the deforming potential and the tidal bulge: 
e{a) = 25{a) = a At. Moreover, for an incompressible body, 
a reasonable assumption for giant planets, the tidal dissi- 
pation function is given by (Goldreich 1963 Efroimsky & 
Williams|[2009l ): 



In order to use the constant time lag model, we must 
consider many values for At. To constrai n this parame- 
ter, we follow the analysis of [Gold reich & Soter ( 1966 ) 
and use the lo-Jupiter system to infer an upper limit for 
Greenberg fc2,p x At^ in giant extrasolar planets. Since Jupiter is 
rapidly rotating, with wj > n (hereafter, J indices refer to 
the value for Jupiter), where n is the orbital mean motion 
of Jupiter's satellites, tidal transfer of angular momentum 
drives the satellites of Jupiter outwards, into expanding or- 
bits. Therefore, the presence of lo in a close orbit provides 
an upper limit for the time lag in Jupiter. Indeed, if Atp 
was too large, the backward evolution of the satellites or- 
bits would imply their disappearance within less time than 
the age of the Solar system, i.e. of Jupiter. For coplanar and 
circular orbits, a dimensionless version of Eq. ([2]) reads: 



Acyc\eE{a) 



tan e{a) 



27r£;peak(^T) 1 - ( f - e(a) ) tan e((7) ' 



(17) 



where Acycic£'(o') is the energy dissipated by the body at 
the frequency a during one tidal period and £'poak(o') the 
maximum energy stored in the perturbation. In the cases 
of interest in the present study, where Q 3> 1, and for non- 
synchronized circular orbits, semi-diurnal tides dominate 
and one can equ al th e average tidal quality factor to the 
one given by Eq. (17) for the relevant frequency 



Q' 



w -k2At\uj-n\. 



(18) 



This formula can be used to estimate the quality factor in 
the case of the jovian planets as long as semi-diurnal tides 
dominate. As the planet tends toward synchronization, the 
dissipative effects of the semi diurnal tides {a — 2\LUp — n\) 
vanish with their frequency. Then, the most dissipative 
tides are the eccentric annual tides {a — n) and 



-k2Atn. 



(19) 



Apart from these two limit cases, no tidal frequency domi- 
nates, and the dissipation is the response of the body to the 
rich spectrum of exciting tidal frequencies. Thus no simple 
relation exists between Q' and At in the general case. 
Although it is tempting to use Eq. ( 19 ) to rewrite the 



tidal equation a nd to keep Q' constant instead of At as done 
by, for exa mple, iMardling fc Lin (|2002 1 , Dobbs-Dixon et al. 
(20041 andlBarker L gilvie| ( |2009[ ), one must keep in mind 
that this procedure is not equivalent either to the constant 
phase lag (i.e. constant Q) or time lag model. Indeed the 
frequency dependence of the phase lag is given by e{a) = 
a I (nQ) and is still proportional to the tidal frequency over 
an orbit as in the constant time lag model, but with a slope 
that is changing during the evolution. 

In S|4] and [Js] we compare the constant time lag model 
with the constant Q' model used by various authors. In or- 
der to allow a direct and immediate comparison with these 
studies, we will choose the values of the couple {Atp, At^,) 



from the relations (fc2Aip 



,k2At^ = 



where (Qp, Q'^,) are the constant normalized quality factors 



used by Miller et al. ( 2009p . This ensures that the effective 
tidal dissipation function is the same in both calculations 
for a given planet with its measured orbital parameters. 



1 



1 



n 



(20) 



where a = a/ao, indices refer to initial values and time is 
counted in units of 



1 



M2 



6 Ms (Ms + Mp)'Rp' GMpfc2,p Aip ' 

where Ms is the satellite mass, t is the typical timescale 
of tidal evolution of the semi-major axis. Injecting angular 
momentum conservation 



Jtnt — 



-a n + GiWi 



Mp + Ms 



into Eq. ( 20 1 , and integrating over time yields 



(21) 



where /3 = (MsMpa§)/((Mp -I- Ms)Cp). Note that this in- 
tegral cannot be performed down to a = because the 
satellite first crosses the corotation radius where the inte- 
grand tends to infin ity, which is an unstable equilibrium 
state for the system ( Hut|1980 ). Furthermore, this result is 
not limited to the case of a planet-satellite system and can 
be directly used to compute the inspiral time of a close-in 
exoplanet (setting p — ★ and S — >■ p) . 

For the lo- Jupiter system, taking t — —4.5 x 10^ yr and 



a{t) equal to the Roche limit in Eq. (21 1 yields kjAtj < 



5 X 10 s. Therefore, for the actual lo- Jupiter system, 

10^, slightly smaller than the 



Eq. (|18|) implies Q\ > 1 x 



valueaerived by Goldreich & Soter ( 1966 1. As discussed by 
these authors, our upper limit on Atp must be multiplied 
by a factor 5 to 7.5, as lo might have been trapped in a low 
order commensurability with Europa and Ganymede dur- 
ing part of its evolution, slowing down the expansion of its 
orbit. This roughly yields 

fcjAtj < 2 - 3 X 10"^ s. 

For sake of easy comparison, we will refer to the quan- 
tity Qq, which is the reduced quality factor computed for a 
reference period of 1 day: 



3Q(27r/lday) 
2 h, 



3 1 day 
22TThhAi' 



The above calculated constraint reads Qq p 1 x 10^. In 
the present calculations, we will examine two cases for the 
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planet under consideration, namely Qq ^ = 10^ and Qq ^ = 
10"^ (A:2,pAtp 2 X 10~2 - 2 x 10"^)^ ^^^i^ taking Q'^[^ in 
the range 10^ - 10'^ (fcz.^At^ ~ 2 x 10"^ - 2 x 10" ^), a 
typical value for solar- type stars ( Ogilvie &: Lin|[2007 |. 

It is important to stress that, if At, or its counterpart 
Q, is poorly known for both planets and stars, its vari- 
ability from one object or configuration to another is even 
more uncertain. For instance, tidal dissipation in planets 
probably differs significantly from the one in brown dwarfs 
because of the presence of a dense co re able to excite iner- 
tial w aves in the convective envelope ( Goodman & Lackner 
20091. Given the highly non- linear behavior of tidal dissi- 



pation mechanisms, the effective tidal dissipation function 
varies not only with the structure of the object or with the 
tidal frequency but also with the amplitude of the tidal 
potential. For example, Q'^ values inf erred from the circu- 
lariza tion of close FGK binary stars (Meibom & Mathieu 
2005 1 , may be lower than the actual Q'^ encountered in 



star-planet systems (Ogilvie & Lin 2007). Consequently, the 
range of values considered here for both Q'^ and Q'p should 
be seen as mean values and be re-evaluated when consid- 
ering specific and/or atypical systems (XO-3, IIAT-P-2 or 
CoRoT-Exo-3 for example). 

4. Effect of the truncation of the tidal equations to 
2'^'^ order in e: Analytical analysis. 



Following the initial studies of Jackson et al. ( 2008 ) , all the 
studies exploring the effect of tidal heating on the internal 
evolu tion of "hot jupiters" (Miller et al. 2009 Ibgui et al. 
2009| have been using a tidal model assummg a constant 
Q value during the evolution. Moreover, in all these calcu- 
lations, the tidal evolution equations are truncated at the 
second order in eccentricity (hereafter referred to as the " 
model"), even when considering tidal evolution sequences 
with non-negligible values of e at earlier stages of evolution. 
Although such a e^-truncated mo del is justified for planets 
and satellit es in the solar system ( Kaula|1963 



Goldreich & 



Soterl [19661), it becomes invalid, and thus yields incorrect 
results for a{t), e{t) and i^tidos for finite eccentricity values. 
The main argument claimed for using such a simple tidal 
model is the large uncertainty on the tidal dissipation pro- 
ces ses in astrophysic al objects. In particular, as detailed 
by Greenberg (2009), the linearity of the response to the 
tidal forcing based on the viscoelastic model may not hold 
in a real object for the large spectrum of exciting frequen- 
cies encountered when computing high order terms in the 
eccentricity. Although the large uncertainty in the dissipa- 
tive processes certainly precludes an exact determination of 
the tidal evolution, it can by no means justify calculations 
which are neglecting dominant terms at finite e. 

Indeed, from a dimensional point of view and prior to 
any particular tidal model, the strong impact of high order 
terms in the eccentricity is simply due to the fact that the 
tidal torque (N) is proportional to {oj — 9)/ {9 being the 
true anomaly) and that over a keplerian orbit, the average 
work done by the torque is of the form 



< N-9 >cx< 



>= 



14e2 



105 4 
4 ^ 



35 6 , _35_ 8 
4 ^ ^128^ 



(1 



o2U5/2 



which is a rapidly increasing function of e (see Appendix 
[a] for the details of the calculation). This means that, al- 
though the mean distance between the planet and the star 



increases with e, the distance at the periapsis strongly de- 
creases, and most of the work due to the tidal forces occurs 
at this point of the orbit. One can see that for e > 0.32 
the high order terms dominate the constant and terms. 
This is a physical evidence that shows that for moderate 
to large eccentricity, most of the tidal effects are contained 
in the high order terms that can therefore not be neglected 
independently of any tidal model. 

In this section we quantify more comprehensively this 
statement. We will demonstrate analytically that: 

— in the context of the Hut model, a truncation of the 
tidal equations at the order can lead not only to 
quantitatively wrong but to quantitatively wrong tidal 
evolution histories, with sequences drastically differing 
from the ones obtained with the complete solution. 

— the rate of tidal dissipation can be severely underesti- 
mated by the quasi circular approximation (e <C 1). 

Furthermore, Q-constant models consider only low order 
terms in obliquity (e;), and thus cannot address the prob- 
lem of obliquity tides and energy dissipation produced by 
this mechani sm. A detailed discus sion on this subject is 



(2009) and will not be reproduced here. 



presented in Levrard et al. (^2007) and Barker & Ogilvie 



4.1. Expanding vs shrinking orbits 

On one hand, considering Eq. ([2]) (with e = for simplifi- 
cation) we can see that, for ujijn < Na{e)/N{e), the tides 
raised on the body i lead to a decrease of the semi-major 
axis, transferring the angular momentum from the orbit to 
the body's internal rotation. It is easy to show that for a 
synchronous planet this condition is always fulfilled, since 
^ ~ ^nXe} ^ ^{e) ^'^^ eccentricity (respectively solid 
and dashed curves of Fig.fTk). As a result, the semi-major 
axis of most short period planets is decreasing. 

On the other hand, truncating Eq. ^ at the order 
for the semi-major axis evolution yields 



1 da 
a dt 



4a 



27 3 ■ 

2 n 



(1 + fe^)^ 
2 n 



(1-1- 23 



(22) 



and the previous condition becomes uji/n ^ (l + 23e^)/(l + 
^e^). Up to 2"^^ order in eccentricity, the pseudo synchro- 
nization angular velocity is given by Weq = (1 -I- 6e^)nFl 
One can see that Weq/n = 1 6 < (1 23 6^)7(1 + ^e^) 
only for 

1 /7 



9 V 2 



0.208. 



This means that even for a moderate eccentricity, e ~ 0.2, 
the truncated model predicts that tides raised on a pseudo 
synchronous planet lead to a growth of the semi-major axis 
instead of a decrease, as obtained by the complete model. 
Therefore, truncating the tidal equations at the order 
for an eccentricity e ^ 0.2 not only predict quantitatively 



^ these equ ations truncated at the or der agree with equa- 
tions in §16 of Ferraz-Mello et al. ( 2008 1, even though they have 



been derived with different methods 
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(b) u/n vs e phase space: truncated equations 



Fig. 1. Pseudo synchronization curve (solid), the d 
curve (dashed) and the e = curve (dotted) for the com- 
plete model (top pannel) and for the truncated one (bottom 
pannel) in the w/nvse phase space. A pseudo synchronized 
planet always follows the w = curve (black solid curve) . It 
always lies in the d < and e < part of the diagram be- 
cause when solving the complete Hut tidal equations, these 
curves do not intersect (top pannel). In contrast, as demon- 
strated in sjij in the 2"'^-order truncated model, the pseudo 
synchronization curve intersects the d = (at e ^ 0.208) 
and e = lines (e ~ 0.326). 



wrong but qualitatively wrong tidal evolutions. The same 
arguments for the evolution of the eccentricity show that 
tides raised on a pseudo synchronous planet lead to a growth 
of the eccentricity for 



e ^ 



0.326 



and not to a decrease. This is illustrated by FigJTJthat shows 
the pseudo synchronization curve (solid), the d = curve 
(dashed) and the e = curve (dotted) for the full model 
(top panel) and the truncated one (bottom panel), in the 
uj/n vs e phase space. As demonstrated before, the pseudo 
synchronization curve crosses the d = and e = lines 
in the 2"^^ order model (Fig.[l|3) whereas it does not when 
solving the complete Hut equations (Fig.jT^). As a result, 
with the truncated model, a pseudo synchronized planet can 



erroneously enter the zone of the phase space where its tides 
act to increase both the semi-major axis and the eccentric- 
ity. While this behavior is not observed with the constant 
phase lag model because it assumes that the star is slowly 
rotating (u-^/n ^1) and that the planet is near synchro- 
nization (ujp/n « 1) - placing them in the d < and e < 
zone of the phase space - this formal demonstration sets 
clear limits on the domain of validity of the quasi circular 
approximation. 

4.2. Underestimating tidal heating 

The key quantity arising from the coupling between the or- 
bital evolution and the internal cooling history of a planet is 
the amount of energy dissipated by the tides in the planet's 
interior, which may compensate or even dominate its en- 
ergy losses. As a result, tides ra ised in an eccentric planet 



can slow down its contraction (Bodenheimer et al. 



Leconte et al. 2009 Baraffe et aT] 20101) or even lead to a 



2001 



transitory phase of expansion ( Miller et al. 2009 Ibgui et al 



2009) . Correctly determining the tidal heating rate is thus 
a major issue in the ev olution of short-period planets. The 



often used formula is (Kaula 1963 
Jackson et al]|2008| ) 



Peale & Cassen 1978 



(23) 



(the 21fc2,p/2(5 is rewritten 63/4g' in Miller et al. (12009)) 



As already stated by Wisdom^ (2008] ) , although this formula 
gives a fair approximation of the tidal dissipation rate for 
the small eccentricity cases, typical in the solar system, it 
severely underestimates the tidal heating for moderate and 
high eccentricities. Fig.[2]illustrates the power dissipated in 
a pseudo synchronized planet as a function of the eccentric- 
it y. It shows that for e > 0.45, the truncat ed formula used 



Miller et al. (2009 ) and Ibgui et al. (2009 ) underestimates 



the actual tidal dissipation rate by more than one order of 
magnitude and by more than a factor 10"^ for e > 0.7, an ec- 
centricity value often advocated by these authors to explain 
the highly inflated planets (see pi). 

From a mathematical point ofview, the fact that a trun- 
cation to 2"'^ order in eccentricity yields such discrepancies 
is due to the presence of (1 — e^)~^^/^ factors in the equa- 
tions for the tidal dissipation. As already stated by [Wisdom 



(2008), for moderate to large eccentricity, this function is 
poorly represented by the first terms of its polynomial rep- 
resentation. Indeed, the first terms of the energy dissipation 
rate are given by: 



tide 



7ifpe2 



1133 . 31845 



28 



224 



381909 
896 



'+0 (el") 



The dissipation rate calculated up to e^" is plotted in Fig 
(dotted curve), where it can be compared with the exact 
result. It is clear that, for e > 0.4, the polynomial devel- 
opments of the tidal evolution equations must be done to a 
much higher degree that done in previous studies, or com- 
plete calculations such as the ones done in Hut\ \198D^ must 
be used. The same argument holds for the evolution of the 
semi-major axis and the eccentricity. Since Eqs. ([2| and ^ 
also contain (1 — e^)~P/^ factors, the decrease of a and 
e is severely underestimated at even moderately large ec- 
centricity when using a 2"*^ order truncated expansion in 
eccentricity. 
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Fig. 2. Tidal energy dissipation rate in a pseudo synchro- 
nized planet (in Watt) as a function of the eccentricity 
calculated with Eq. (13 1 {solid curve) and with the trun- 
cated formula (Eq. (23l; dashed). The ratio of the two 



curves only depends on the eccentricity and not on the 
syste m's parameters. For e = 0.45, the approximation 
(Eq. (23)) underestimates the tidal heating by a factor 10. 
The actual values were derived using HD 209 458 b param- 
eters: Mp = 0.657Mj, flp = 1.3 2i?j, Ah = 1.101 M©, 



a = 0.047 AU ( [Knutson et al.|2007j ). Q' = 10^ (see The 



dotte d cu rve gives the dissipation rate calculated up to e^'^ 
(Eq.m). 



5.1. Calculations at low eccentricity 
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In particular, as discussed in the next section, a high ec- 
centricity (e > 0.6) cannot be maintained for a few 100 Myr 
to a few Gyr i n a system like HD 2 09 458 in agreement with 

se e Fig . [4] below). This 
( 2009| who find that 



the results of Miller et a l 
is in contrast with |Ibgui 



2009) 
& Burrows 



the radius HD 209 458 b can be matched and that the sys 
tem can sustain a significant eccentricity up to the observed 
epoch. Such discrepancies between these two studies based 
on the same tidal model may reveal differences in the im- 
plementations of the tidal equations, or a difference in the 
calculation of interior structures or boundary conditions. 



5. Effect of the truncation to 2'^'^ order in e: 
Simulation results 

In this section, we present the comparison of the results of 
our complete model with the "e^ model". We have calcu- 
lated evolutionary tracks of the tidal evolution for various 
transiting systems, coupling the internal evolution of the 
object either with our tid al m odel or with the "e^ mo del" 



used in [Miller et al] (|2009) and |Ibgui fc Burrowsj ( [20091 ). In 
order to ensure a consistent comparison with these authors, 
we directly convert their set of tidal parameters. Since our 
model assumes a constant time lag, and not a constant Q' 
value, a history track computed with the Q' "e^ model" 
with a constant couple (Qp, Q'i,) is compared with a his- 
tory track computed in our model with a constant cou- 
ple (fc2Atp, /c2Ai*) given by {k2Atp = 2„2sQ' ' ^sAt^ = 

Q, ) (See ^and Eq. ( 19 )). This ensures that - although 
our calculations are conducted with a constant At - the 



quality factor computed with Eq. ( 19 ) in the object at the 



present time is the same as the one used in the Q constant 
model. 



Fig. 3. Consistent tidal/thermal evolution of TrES-l b 
computed with our constant time lag model (solid line) 
and with the "e^ model" (dashed line). This is a 0.7 6 Mj 



planet orbiting a 0.89 Mq star ( [Winn et aL][2007[ ) 
error bars are the measured paranieters with the la uncer 



tainty . The dashed curve is comparable to Fig. 7 of [Miller] 
et al. (2009) and was computed with the same parameters 
(g'p = 10°-^, Q'^ = 10^). As expected, in the low e limit, 



the tidal dissipation rate is well approximated by Eq. ( 23 ) 
and the two models yield similar evolutions, although the 
merging time depends on the rheology used. 



We first compare the results of the two models on a 
system which has a zero measured eccentricity and is not 
inflated, namely TrES-l. Such a system does not require a 
substantial initial eccentricity for its observed properties to 
be reproduced and thus provides an opportunity to t est the 
quasi c ircula r limit, where the "e^ model" used by Miller 
al ' ■ 



^ model" (dashed 
As expected, in 



et al.[ fl2009| and our model should yield similar results 
Fig. [3 illustrates the results of the integration of the coupled 
internal/orbital evolution equations with our constant time 
lag model (solid curve) and with the "e 
curve) for an initial eccentricity of 0.07, 
this low eccentricity limit both models yield very similar 
tracks: the eccentricity is damped to zero in a few Gyr and 
the semi-major axis decreases until the planet reaches the 
Roche limit and merges with the star, because the system 
does not hav e enough angular momentum to reach a stable 



2009). In this case, 



equilibrium ([Hut[|198(M Levrard et all 
tidal heating is not sufficient to significantly affect the ra- 
dius of the planet which keeps shrinking steadily as it cools. 
Note, however, that, although the qualitative behavior of 
the evolution is the same, the hypothesis made on the rheol- 
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ogy of the body can influence the age at which the merging 
occurs. 



5.2. Calculations at high eccentricity 
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Fig. 4. Consistent tidal/thermal evolution of X0-4b {thin, 
black) and HD 209 458b {thick, blue) computed with our 
constant time lag model {solid line) and with the "e^ 
model" {dashed lin e). X0-4b is a 1.72 Mj planet orbit- 
ing a 1.32 Mq star ( [McCullough et al.|2008| ). HD 209 458 b 
is a 0.657 Mj planet orbiting a 1.01 Mq star ( Knutson 



et al. 2007 1. The dashed cu rves are comparable to Fig. 8 
and 10 of Miller et al. ( 2009 ) and were computed with the 
same parameters {Q'p = 10 , Q'^ — 10^). In this large ec- 
centricity regime, using the same quality factor, the "e^ 
model" underestimates the tidal dissipation rate by 2 or- 
ders of magnitude and thus overestimates the star-planet 
merging timescale by a factor 10 to 10'^. 



In the moderately to highly eccentric regime, the tidal 
dissi patio n rate can no l ong er be approximated by Eq. ( 23 ) 
(see ^4.2). Instead, Eq. (13) must be used and yields - as 
shown by Fig. [2]- a much more important dissipation rate. 
As a result, tidal evolution takes place on a much shorter 
time scale, and both the eccentricity damping and the merg- 
ing with the star occur earlier in the evolution of the planet. 
For illustration. Fig. |^ portrays the possible thermal/tidal 
evolution (for given initial conditions) for X0-4b (thin 
black curves) and HD 209 458 b (thick blue curves) com- 
puted with the "e^ model" (dashed) and with our model 
(solid) . The dashe d curves are s imilar to the ones displayed 
in Fig. 8 and 10 of Miller et al. ( 2009 1 . As mentioned above 
and illustrated in Fig.l4|l, the energy dissipation is much 



larger when fully accounting for the large eccentricity. The 
evolution of the planet can exhibit two different general 
behaviors: 

— The planet first undergoes a phase of contraction and 
rapid cooling before the tidal heating due to the large 
initial eccentricity starts to dominate the energy bal- 
ance of the object, leading to a phase of radius infla- 
tion (as shown by Fig. [5] for a test case). This speeds 



o 



-10 




10- 



0.01 



Age (Gyr) 

Fig. 5. Internal energy balance in the evolving planet. Solid 
line: luminosity of the object with tidal heating. Dotted 
line: luminosity of the object without tidal heating. Dashed 
line: tidal energy dissipation rate. The object contracts as it 
cools until the energy input balances its thermal losses and 
sustains a higher entropy in the gaseous envelop, yielding 
a larger radius. 



up the damping of the eccentricity and the decrease of 
the semi-major axis, since a and e oc R^. When the 
eccentricity becomes small enough, a "standard" con- 
traction phase begins and lasts until the planet mer ges 
with the star (due to stellar tides; Levrard et al. 2009 1 or 

- if enough angular momentum is present in the system 

- until both tidal and thermal equilibria are achieved. 
Such a beha vior has already bee n identified by Miller] 
et al.| (|2009[) and [Ibgui fc Burrows| ( 2009 ) but, because 



these authors use truncated tidal equations, they find 
that a large eccentricity can be maintained for a few 
Gyr and keep inflating the planet at a late time, as il- 
lustrated on Fig. |4] (dashed curves); whereas it is not the 
case. 

— In some extreme cases, such as HD 209 458, for the 
initial conditions corresponding to the ones in Fig. 4] 
the tidal heating can overwhelm the cooling rate of the 
planet by orders of magnitude and lead to a spectacular 
inflation of the planet and thus to a rapid merging with 
the star. This stems from a combination of different ef- 
fects. First of all, as mentioned above, the expansion of 
the radius accelerates the tidal evolution and thus the 
decrease of the orbital distance. Furthermore, the Roche 
limit (flR = aRp -y/M^/Mp where a is a constant which 
depends on the structure of the body and is equal to 
2.422 for fluid objects) increases with the radius of the 
planet, extending the merging zone. 

As clearly illustrated by these calculations, using tidal 
equations truncated at second (e^) order leads to severely 
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erroneous evolutionary tracks for initially moderately (e > 
0.2) or highly eccentric systems. Indeed, the complete tidal 
model shows that, for the initial conditions and Q pa- 
rameter values cho sen by Ibgui & Burrows (20091 and 
Miller et alT] ( |2009[ ), HD 209458 b would in factTiave dis- 
app eared! As mentioned earlier, the use of such a quasi 
circular approximation cannot be justified by the uncer- 
tainty on the quality factor, as the discrepancy in the char- 
acteristic evolution timescales can amount to 3 orders of 
magnitude in some cases, depending on the initial eccen- 
tricity. Conversely, trying to infer values for the stellar or 
planetary tidal quality factors Q from tidal evolution cal- 
culations performed with the truncated model will lead 
to severely inaccurate values. 

6. Global view of transiting systems 

As mentioned earlier, tidal heating has been suggested by 
several authors to explain the anomalously large radius of 
some giant close-in observed exoplanets. As demonstrated 
in ^|5] the previous calculations, all based on constant-Q 
models truncated at the order yield inaccurate results 
when applied to significantly (initial or actual) eccentric 
orbits - a common situation among detected exoplanetary 
systems. In this section, we revisit the viability of such a 
tidal heating mechanism to explain the large observed Hot 
Jupiters radii with the present complete Hut tidal model. 
We first examine the properties of the known transiting 
systems. Then we show that, although indeed providing a 
possible explanation for some transiting systems, the tidal 
heating hypothesis fails to explain the radii of extremely 
bloated planets such as - among others - HD 209 458 b, 
TrES-4b, WASP-4b or WASP-12b, in contrast with some 
previously published results based on truncated tidal mod- 
els (See 0. 

It is now well established that a large number of transit- 
ing giant exoplanets are more infiated than predicted by the 
sta ndard cooling theory of irradiated ga seous giant planets 



(see Udry fc Santos|2007 



Baraffe et al. 2010 for reviews). In 



order to quantify this effect, we have computed the radius 



predicted by our standard model, described in §2.1[ for the 
54 transiting planets detected at the time of the writing 
of this paper, with Mp > O.SAfj (about a Saturn mass). 
We define the radius excess as the difference between the 
observed radius and the one predicted by the model at the 
estimated age of the system, denominated i?irrad- Results 
are summarized in Fig.|6] The existence of objects below 
the R = i?irrad hue is a clear signature of the presence 
of a dense core and/or of th e enrichment of the gaseous 
envelop (Baraffe et al. 2006; Fortney et al. 2007; Burrows 



et al.||2007; ,Baraffe et al.„20087]reconte et al.| |2009). Note 
that most of the objects significantly below tliis line are 
in the M < IMj region, and can be exp lained with a 
Mz/Mp>0.10 heavy material enrichment (Baraffe et al. 



20081, in good agreement with predic tions of the core ac 
cretion scenario for planet formation ('Baraffe et al 



2006 



Mordasini et aL|[2009 ). Interestingly enough, all the planet 
radii in the R ^ -Rirrad region of Fig.|6] show no significant 
eccentricity and can be explained by including the presence 
of a core in their internal structure and an orbital evolution 
with a low initial eccentricity, independently of the chosen 
tidal parameters. 

Among the 39 remaining objects, we will focus on the 
most extremely inflated ones to investigate the validity of 



5- 








0.2 
















1 


0.0 






-0.2 






Fig. 6. Relative radius excess between the observationally 
and the theoretically determined values for 54 transiting 
systems. Objects above the dashed line are considered as 
being anomalously bloated compared with the prediction 
of the regular evolution of an irradiated gaseous planet. All 
the objects below this line can be explaine d by a heavy 
material enrichment in the planet's interior (Baraffe et al. 
20081). 



the tidal heating hypothesis to explain their abnormally low 
density, as they provide the most stringent cases to examine 
the viability of this scenario. For sake of simplicity and 
to avoid introducing further free parameters in our tidal 
model, and since our aim is to derive an upper limit for the 
radius that a planet can achieve under the effect of tidal 
heating, we will not consider the presence of heavy element 
enrichment in our calculations. 

Our calculations proceed as follows: 

1. For each of the systems, a range of initial semi-major 
axis ([ai,min, ai,max]) IS fouud by a backward integration 
of the tidal equations, from present-day observed values. 

2. Evolutionary tracks, coupling consistently the 
gravothermal evolution of the irradiated planet 



and the tidal heating source (Eq. (13)), are then com- 
puted for fli S [ai_min, fti.max] and an initial eccentricity 
Ci £ [0,0.8]. The plausibility of such initial conditions 
as a remnant of early plane t-disk and/or planet-planet 



interaction is discussed in Miller et al. (20091. Since 



3. 



total angular momentum is conserved during the tidal 
evolution, the initial spin rate of the star is calculated 
by satisfying the equality between the initial and 
the presently observed value of the system's total 
angular momentum. Calculations are performed with 
Q'g^^ = 105 and 10^ and Q'^ ^ = 10^ and 10^ (See ^ioi 
a detailed discussion). 

For each evolutionary calculation, the departure 
from a given measured quantity is defined as 
Sx:{t) = ^ ^i^y-^p ^ ^ where x refers to a, e, Rp, or 

uJif and ax to their measured uncertainty. When no er- 
ror bar has been measured for the eccentricity and e = 
has been assumed in the light curve analysis, we take 
(Te — 0.05. We consider that the evolution accurately re- 
produces the presently measured data if there is a time 
interval (compatible with the age of the system) within 
which all the Sx^s are smaller than 1, meaning that each 
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one of these parameters agrees with the measured one 
within 1 a. 
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Fig. 7. Set of initial conditions yielding evolutions consis- 
tent with the actual orbital parameters of HD 209458 b. 
These runs assume Qq p = 10^ and Qg^ = 10^. As large 
eccentricity speeds up the tidal evolution, when initial ec- 
centricity of the orbit increases, the initial semi-major axis 
must also increase to recover the observed parameters at 
the age of the system. 



Fig.[7] portrays a grid of evolution history initial con- 
ditions that are found to be consistent with the observed 
parameters of HD 209 458, at the age of the system. As 
expected, an initially more eccentric system must have a 
greater initial separation to end up at the same location. 
This stems from the fact that \a\ is a monotonically in- 
creasing function of e for a slowly rotating star (as obtained 
from Eq. ([2| for Wp = Wcqu and w^/ri ^ 1). Although, for 
these extremely bloated planets, we do find evolutionary 
tracks that lead to the presently observed orbital param- 
eters, none of these solutions can reproduce the presently 
observed radii, as illustrated on Fig. [8] for the case of HD 
209 458 b. Indeed, the major phase of ecce ntricity damping, 

occurs too early 
the tidal heating 



as given by Eq. Q and discussed in j ]5.2 
in the evolution, so that a large fraction o: 
energy dissipated in the planet has been radiated away at 
the age of the system (~ a few Gyrs) and can no longer 
provide enough energy to slow down gravitational contrac- 
tion. The same behavior is found for other bloated systems 
such as WASP-12, TrES-4 and WASP-4 whose best evo- 
lutionary tracks are shown in Fig.|9] For all these systems, 
although a solution matching the presently observed orbital 
parameters can be found, tidal dissipation occurs too early 
to reproduce the present values of the planet radii. 

These results, based on complete tidal evolution calcu- 
lations, show that the tidal energy dissipated in the planet's 
tidal bulges, although providing a viable explanation to the 
large radius of many short-period planets (such as OGLE- 
TR-211b shown in Fig. 10), is not sufficient to explain the 
radii of the most bloatea planets, at the age inferred for 
these systems. In that case, an extra mechanism, besides 
tidal heating, must be invoked to solve this puzzling prob- 
lem. Surface winds driven by th e powerful incident stellar 
flux ( Showman fc Guillot||2002 1 , converting kinetic energy 
to heat by dissipation within the tidal bulge and thus reach- 
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Fig. 8. Consistent tidal/thermal evolution of HD 209458 b 
with different initial conditions (solid and dashed) com- 
puted with our constant time lag model. HD 209 458 b is 
a 0.6 57 Mj planet orbiting a 1.01 Mq star ( Knutson et al. 
2007[). These runs assume Q'g p = 10^ and Q'^ ^ = 10°. For 



comparison, the radius and luminosity of an isolated planet 
(no tidal heating) is shown on the lower panels (dotted 
curves) . Even though these evolutions recover the presently 
observed orbital parameters for the system, the eccentricity 
damping arises too early during the evolution, leading to 
insufficient tidal dissipation at present epoch to explain the 
inflated radius. 

ing deep enough layers to affect the planet's inner isen- 
trope, or inefficie nt large-scale convection d ue to a com- 
position gradient (Chabrier &: Baraffe 2007) could be the 
other mechanisms to be considered with tidal dissipation 
to lead even tually to these large planet radii (see Baraffe 
et al. (2010) for discussion). 



7. Discussion and conclusion 

In this paper, we have demonstrated that the quasi-circular 
approximation (e ^ 1, i.e. tidal equations truncated at the 
order e^) usually made in tidal calculations of transiting 
planet systems and valid for our Solar system planets, is 
not valid for the exoplanetary systems that have - or were 
born with - an even modestly large (e>0.2) eccentricity. 
As shown in f|4j although the real frequency dependence 
of the tidal effect remains uncertain, there are dimensional 
evidences that for eccentric orbits, most of the tidal effect 
is contained in the high order terms and that truncating 
the tidal equations at 2"*^ order in eccentricity can overes- 
timate the characteristic timescales of the various orbital 
parameters by up to three orders of magnitude. Therefore, 
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Fig. 9. Evo lutionary tracks for WASP-12b (solid, Hebb 



et al | |2009[ ), TrES- 4b (dashed , Dacm gen et"al||2009p aUd 
WASP-4b (dotted, [Winn et al.|2009 a) that lead to the best 
agreement with the observed orbital parameters for these 
systems. These runs assume Q'q „ = 10^ and Q;,^ = 10^. 
Tidal dissipation is not sufficient to sustain the large radii 
observed for these planets. 



Fig. 10. Evolutionary tracks for OGLE-TR-211 b that lead 
to the best agreement with the observed p arameters, is a 
L03 Mj planet orbiting a 1.33 Mq star (Udalski et al. 
2008[ ). These runs assume Qf, p = 10^ and Q'^^ = 10°. The 
radius evolution of OGLE-TR-211 b without tidal heating is 
plotted in panel (c) for comparison (dashed). For this mod- 
erately inflated planet, tidal heating is sufficient to sustain 
a large radius for the age of the system. 



truncating the tidal equations at the second order can by 
no means be justified by invoking the large uncertainty in 
the dissipative processes and their frequency dependence. 
Therefore, high order tidal equations should be solved to 
derive reliable results for most of the existing exoplanet 
transiting systems. This need to solve the complete equa- 
tions is met by any tidal model. In this context, even though 
no tidal model can claim describing perfectly a two body 
evolution, we recall that the Hut model is at least exact in 
the weak friction viscous approximation (see Sjs]). 

We have tested our complete tidal model on several in- 
flated planets to find out whether or not tidal heating can 
explain the large radius of most of the observed transit- 
ing systems. Although this mechanism is indeed found to 
be sufficient to explain m ode rately bloated planets such as 
OGLE-TR-211 b (see Fig.[lO]), we have been unable to find 
evolutionary paths that reproduce both the measured ra- 
dius and the orbital parameters of HD 209 458 b, WASP- 
12 b, TrES-4b, and WASP-4b (see Fig.js] and Fig.]! for 
their inferred age range. The main reason is the early cir- 
cularization of the orbit of these systems. As demonstrated 
in the paper, this stems from the non-polynomial terms in 
eccentricity in the complete tidal equations, which are miss- 
ing when truncating the equations at small e-order. The 
present results, based on complete tidal equations, show 
that tidal heating, although providing an important con- 
tribution to the planet's internal heat budget during the 



evolution, cannot explain alone the observed properties of 
all exoplancts. 

This is in contrast with some of the conclusions reached 
in previous studie s. Ba sed on truncated t idal models, |lbgui| 
'& Burrows (2009) and Ibgui et al. (2009) find evolutionary 
tracks that match observed parameters for HD 209 458 b, 
WASP-12b, and WASP-4b and thus suggest that the tidal 
heating is the principal cause of the large radii of Hot 
Jupiters. 

These particular properties of Hot Jupiters, including 
the extreme cases of the most severely bloated planets, can 
only be explained if the following explanations/mechanisms 
occur during the system lifetimes: 

— Early spin up of the star: si mulations of the rotational 
evolution of solar-like stars (Bouvier et al. 1997) show 



that after the dispersion of the accretion disk, the rota- 
tion rate of the contracting star increases due to angular 
momentum conservation, until magnetic braking takes 
over. Considering Eq. we see that stellar tides act 
as an eccentricity source if ^ ^ y| ffj^ ■ Investigating 
whether the duration of this phase lasts long enough 
and whether the magnitude of this effect is large enough 
to drive enough eccentricity requires performing consis- 
tent star/planet thermal/tidal calculations and will be 
investigated in a forthcoming paper. 
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Presen ce of a third body: as proposed by |Mardling| 



(2007), a low mass terrestrial planet can drive the ec- 
centricity of a massive giant planet during up to Gyr 
timescales. Accurate enough observations are necessary 
to support or exclude the presence of such low-mass 
companions. 

— As mentioned earlier, combining tidal heat dissipation 
with other mechanisms such as surface winds, due to the 
stellar insolation, dissipating deep enough in the tidal 
bulges, or layered convection within the planet's interior 
may provide the various pieces necessary to completely 
solve the puzzle. 

In conclusion, the suggestion that tidal heating is the 
main mechanism responsible to solve the problem of anoma- 
lously large short-period planets, as sometimes claimed 
in the literature, must be reformulated more rigorously: 
although providing a non-negligible contribution to hot- 
Jupiter heat content, tidal dissipation does not appear to 
provide the whole explanation. Further studies are thus nec- 
essary to eventually nail down this puzzling issue. 
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Appendix A: Tidal evolution equations for finite 
eccentricity and obliquity. 

The present calculation of the tidal ev olution equations ex- 
tends the formulas given in |Hut ( 1981 1 to any obliquity. We 
consider a system of two deforrnable bodies of mass Mi and 
M2. The demonstration follows three main steps. First we 
compute a vector expression for the tidal force and torque. 
Second we derive the variation of the rotation rate, obliq- 
uity and orbital angular momentum thanks to this expres- 
sion of the torque and using the total angular momentum 
conservation. Finally the evolution of the semi-major axis 
and eccentricity are obtained from the expression of the en- 
ergy dissipated by tides in the dcformable body. The total 
amount of energy dissipated by tides in one of the bodies 
is a direct product of the calculation. 

Up to the quadrupolar terms in the tidal deformation, 
the mutual interaction of the tidal bulges is negligible and 
we can consider separately the effects of the tides raised in 
each body and sum them up at the end of the calculation. 
Let us consider the effect of the tides raised in a dcformable 
body (say Mi , hereafter the primary) in interaction with a 
point mass (say M2 the secondary). The mass distribution 
of a dcformable body in a quadrupolar tidal potential can 
be mimicked by a central mass Mi — 2m and two point 
masses at the location of the tidal bulges (r+, r_) of mass 
m with ||r+|| = ||r_|| — Ri the radius of the primary and 

m = \k2M2 ("tt) where fc2 is the love number of degree 2 
of the primary and r is the distance between the center of 
the two objects. Since we consider a constant time lag Aii 
between the deforming potential and the tidal deformation 
in the frame rotating with the primary, f^- — r(t — Ati) ( 
f refers to the unit vector associated to r) in this frame. 
Let 9 be the orbital rotation vector collinear to the orbital 
angular momentum and whose value is the instantaneous 
variation rate of the true anomaly 9 of the bodies in their 
keplerian motion and wi the rotation vector of the primary. 
Thus, to first order in Ati, 



f+ = f (< — Ail) ~ f — A^i f X (wi 



(A.l) 



the amplitude of the tidal bulges also lags behind the de- 
forming potential and is given by 

1 / R " 

m(t) = -k2M2 1 —, \ — r 

^ ' 2 \r{t- Ml) 

^Mfs (^^^(l + S^^Ati), (A.2) 

and the force exerted by this mass distribution on the sec- 
ondary is 

GM1M2 Gk2MlR'^ , r , ^ ^ 

F = ^i-^-f - 3 — ^..^ n i + 3-Aii}-f 



ak2MlR^, ^ , ^ 
3 — ^ AtiiuJi -9) xr. 
r' 



(A.3) 
(A.4) 

and the angular momentum conservation yields 

N = h = -L, (A.5) 



Thus the tidal torque reads: 

Mk2MlR\ ^ , ■ ^ 

N = 3 — „^ ^ f X Ati(a;i - 9) x r 



where h — j^^'^^ na^Vl — is the orbital angular mo- 
mentum and L = CiUJi the rotational angular momentum 
of the primary. We can then simply derive the rate of an- 
gular velocity variation: 



dt 



(Ciui) = i = L-L = -N-L 



(A.6) 



This product can be carried out by projecting in any 
base. We choose the base defined by h = (0, 0, h) and 
wi = (wisinei, 0, wicosei) where ei is the obliquity. In 
this base. 



Ck.M'^R^ I wisineisin2(6'-f V) 

N = 3 ^ 7 ^ Ati -uJi sinei cos(6l + V') sin(6l + V) 



LUi COS El — 



(A.7) 



where V' is the longitude of the periapsis in this base. The 
precession of the periapsis occurring on a much shorter 
timescale than the tidal migration, we can average the tidal 
torque over ip. This yields: 



N = 3^^A.i 



^LOi smei 


UJi cos El — ( 



(A. 



We can compute the dot product in Eq. (A.6) giving (with 
xi — cosei) 



dCiwi _ Gk2AtiMlR\ 
dt ^6 1 ^ 



l + xl 



LUl 



(A.9) 



The mean rotation rate variation (Eq. ([9])) is obtained by 
averaging over a keplerian orbit using 



1 + e cos ( 



1 



Torh Joibit ^ ^ 



• dt = 



■bit 



9r^ 



(A.IO) 



d9^Vt{e) (A.ll) 



and 



Torb Jorbit ^ ^ 



dt = nN{e) 



(A.12) 



where dt = d9/9, and the angular momentum conservation 
over one orbit is used to express 9 with respect to 9 (see 
Eq. (Ill]) and Eq. (fsl) for the definition of n{e) and N{e)). 



The variation of the obliquity can be obtained with: 

ii = (L-h)=L-h + L-h. (A.13) 



Carrying out the differentiation and using Eq. (A.5) yields 
(L-h)(N-L) (N-h) 



Xi 



i2 L 

(h-L)(N-h) (N-L) 



(A.14) 



Subsituting Eq. (A.8) for N we get after simplification 
dei _ 3 Gk2AtiMlRl . ' 



dt 



sm£i 



^1 o » 
Ci Giui 



h 



(A.15) 
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Averaging over an orbit using Eqs. ( A.ll |-( A.12 ) gives 
Eq. pO|. 



To obtain the variation of the semi-major axis and ec- 
centricity, we must compute the work done by the tidal 
force on the secondary: 



< Eoih >- 



1 



F • dr, 



1 

Torb 



orb J orbit 

{rFr +N 



Idt, 



(A.16) 



(A.17) 



■bit 



where is the radial component and N • 9 the normal one 



The radial forces in r ^ and r ^ in Eq. (A.3) are conser 



vative and do not dissipate energy. The radial part of the 
work is computed using 



: smt 



and thus 
1 



Toi h J orbit 



■dt = 



(A.18) 



(A.19) 



(see Eq. ([S]) for the definition of Ne{e)). The normal com- 
ponent can be written 



One can see that the dissipa ted energy is positive for any 
value of e and xi as expected (Hut 1981 ), and that it is min- 
imum when the body is pseudo synchronized. S ubst ituting 
i jji by the pseudo synchronization rate (Eq. ( |12[ )), Eq. 
(A. 26 1 simplifies to Eq. (13) that can be used for a close- 



in gas-giant exoplanet. For rocky planets locked in syn- 
chronous rotation by their permanent quadrupolar mass 
distribution, t he hea ting rate can be estimated setting 
wi = n in Eq. (|A.26|). 

In fine, the complete equations taking into account tides 
in both bodies are obtained by computing the effects of the 
tides raised in the secondary (given by the same equations 
with 1 <=^ 2) and by adding them up to the effects of the 
tides in the primary. 



2 r>5 / - e 



d0. 



(A.20) 



Again, averaging is carried out using Eqs. (A.ll )-( A.12 ). 
After integration, 



< Sorb >- 2/^1 



N{e)x: 



Naie) 



(A.21) 



The variation of semi-major axis due to the tides raised in 
the primary (Eq. ([2])) is directly given by 



< Eorh > — 



d GM1M2 GM1M2 



dt 2 a 2 a? 

Since the orbital angular momentum is given by 



h=,lG^^ail-e^). 



Ml + M2 



(A.22) 



(A.23) 



the variation of the eccentricity can be obtained by differ- 
entiating h with respect to t: 



2h 
~h 



2ee 



(A.24) 



Only total angular momentum is conserved, t hen h — 
-d{Cl l^Jl)/dt and substituting Eq.p]), and ( |A.23[ ) in 
(A.24) gives one of the two terms ofTlq. |6| corresponding 
to the contribution of one of the two bodies (star or planet) 
for the evolution of the eccentricity. Finally, the rate of tidal 
energy dissipation into the primary is 



d GM1M2 d . 



(A.25) 



Thus, substituting Eqs. ( |A.9| ) and ( |A.21[ ) in Eq. ( |A.25[ ) 
gives 

-Etidcs — 2/^1 



'N^{e)-2Nie)x,^-^+(^f)m{^)' 



(A.26) 
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